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We consider a billiard model of a self-bound, interact- 
ing three-body system in two spatial dimensions. Numerical 
studies show that the classical dynamics is chaotic. The corre- 
sponding quantum system displays spectral fluctuations that 
exhibit small deviations from random matrix theory predic- 
tions. These can be understood in terms of scarring caused by 
a 1-parameter family of orbits inside the collinear manifold. 



Billiards are interesting and useful models to study the 
quantum mechanics of classically chaotic systems jUI). 
In particular, the study of the completely chaotic Sinai 
billiard || and Bunimovich stadium Q showed that the 
quantum spectra and wave functions of classically chaotic 
system exhibit universal properties (e.g. spectral fluctu- 
ations) H as well as deviations (e.g. scars of periodic or- 
bits) |(| when compared to random matrix theory (RMT) 
predictions of the Gaussian orthogonal ensemble (GOE) 

While the theory of wave function scarring has reached 
a mature state in two dimensional systems p[-pT| there is 
a richer structure in more than two dimensions. In par- 



V(r) = 



for r < a, 
oo for r > a. 



(2) 



ticular, invariant manifolds in billiards 12 and systems 
of identical particles (l3) may lead to an enhancement in 
the amplitude of wave functions jl4|,[ll| provided classical 
motion is not too unstable in their vicinities. 

The purpose of this letter is twofold. First we want 
to quantize a billiard model for three interacting par- 
ticles and study a new type of wave function scarring. 
Second, our numerical results strongly suggest that the 
system under consideration is chaotic and ergodic. This 
is interesting in view of recent efforts to construct higher- 
dimensional chaotic billiards |p^|-|l8|. 

This letter is organized as follows. First we introduce 
a billiard model of an interacting three-body system and 
study its classical dynamics. Second we compute highly 
excited eigenstates of the corresponding quantum system 
and compare the results with RMT predictions. 

Recently, a self-bound many-body system realized as 
a billiard has been studied in the framework of nuclear 
physics We want to consider the corresponding 

three-body system with the Hamiltonian 
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where fi is a two-dimensional position vector of the i-th 
particle and pi is its conjugate momentum. The two- 
body potential is 



The particles thus move freely within a convex billiard in 
six-dimensional configuration space and undergo elastic 
reflections at the walls. Besides the energy E, the total 
momentum P and angular momentum L are conserved 
quantities which leaves us with three degrees of freedom. 
In what follows we consider the case P = 0, L = 0. 
To study the classical dynamics it is convenient to fix 



the velocity vf + u| + = 1 and perform computations 
with the full Hamiltonian (|l|) without transforming to 
the subspace P = 0, L = 0. We want to compute the 
Lyapunov exponents of several trajectories. To this pur- 
pose we draw initial conditions at random and compute 
the tangent map |2(],[l7| while following their time evolu- 
tion. To ensure good statistics and good convergence we 
follow an ensemble of 7 x 10 4 trajectories for 10 5 bounces 
off the boundary. All followed trajectories have positive 
Lyapunov exponents. The ensemble averaged value of 
the maximal Lyapunov exponent and its RMS deviation 
are Xa — 0.3933 ± 0.0015, while the second Lyapunov 
exponent is also always positive. Thus, the system is 
chaotic for practical purposes. However, we have no gen- 
eral proof that no stable orbits exist. The reliability of 
the numerical computation was checked by (i) comparing 
forward with backward evolution, (ii) observing that en- 
ergy, total momentum and angular momentum are con- 
served to high accuracy during the evolution and (iii) 
using an alternative method to determine the Lya- 
punov exponent. 

The considered billiard possesses two low-dimensional 
invariant manifolds that correspond to symmetry planes. 
The first "collinear" manifold is defined by configurations 
where all three particles move on a line. The dynamics 
inside this manifold is governed by the one-dimensional 
analogon of Hamiltonian (|j). After separation of the 
center-of-mass motion one obtains a two-dimensional bil- 
liard with the shape of a regular hexagon. This system 
is known to be pseudo-integrable p2| , p3[ . To study the 
motion in the vicinity of the collinear manifold we com- 
pute the full phase space stability matrix for several pe- 
riodic orbits inside the collinear manifold which come in 
1-dim families and can be systematically enumerated us- 
ing the tiling property of the hexagon. All considered 
types of orbits except two are unstable in the transverse 
direction: (i) The family of bouncing ball orbits (i.e. two 
particles bouncing, the third one at rest in between) is 
marginally stable (parabolic) in full phase space, (ii) The 
family of orbits depicted in Fig. |] is stable (elliptic) in two 
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transversal directions and marginally stable (parabolic) 
in the other 10 directions of 12-dim phase space. Though 
this behavior does not spoil the ergodicity of the billiard 
one may expect that it causes the quantum system to 
display deviations from RMT predictions. Note that this 
family of periodic orbits differs from the bouncing ball 
orbits which have been extensively studied in two- and 
three-dimensional billiards |24|,[l6| since (i) it is restricted 
to a lower dimensional invariant manifold, and (ii) it is 
elliptic (complex unimodular pair of eigenvalues) in one 
conjugate pair of directions. 

The second invariant manifold is defined by those con- 
figurations where two particles are mirror images of each 
other while the third particle is restricted to the motion 
on the (arbitrarily chosen) symmetry line. Inside this 
manifold one finds mixed (i.e. partly regular and partly 
chaotic) dynamics. However, the motion is infinitely un- 
stable in the transverse direction due to non-regularizable 
three-body collisions. 

The quantum mechanics is done using the coordinates 
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Here p, 9' and <p' describe the intrinsic motion of the 
three-body system while x and <j> are the center of mass 
and the global orientation, respectively. In a second 
transformation we apply a rotation of 7r/2 around the ab- 
scissa corresponding to spherical coordinates (p,6',cp'), 
namely t&nip = — cot 8'/ cos <p' and cos 9 = sin 9' simp' 
and obtain for the Laplacian in the subspace P = 0, L = 
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Products of Bessel functions and spherical harmonics 

^k Mu {pAv) = {kp)- 1 J2i+i{kp)Y^(8 lV ) (5) 
are eigenfunctions of the Laplacian (||) 

Aip k ,l,h (p, 0, <p) = -tffaM, (p, 9, tp), (6) 

with the usual relation between wavevector and energy 
k = h^ 1 (2mE) 1 ^ 2 . Fig. p| shows a picture of the billiard 
taking (p, 9, ip) as spherical coordinates. The billiard 
possesses a D^h symmetry. In the fundamental domain 
(0, <p) S (0, 7r/2) x (— 7r/6,7r/6) the boundary is given by 
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The collinear manifold is the equatorial plane 9 — n/2. 
The second invariant manifold is given by the vertical 
symmetry planes <j> = ±tt/6. Note that in this represen- 
tation, classical geodesies of the billiard between two suc- 
cessive collisions are not straight lines since the centrifu- 
gal potential is stronger than in Euclidean case. In what 
follows we restrict ourselves to the fundamental domain 
and choose basis functions that fulfill Dirichlct boundary 
conditions. These are bosonic states. 

We are interested in highly excited eigenstates. These 
may be accurately computed numerically by using the 
scaling method developed in ref. |^5j and applied to a 
three-dimensional billiard by one of the authors [fj4| . 
This method works efficiently only when a suitable posi- 
tive weight function is introduced in a boundary integral. 
To this purpose we note that the radial part of (Q) looks 
like a 4-dim Laplacian. Extending the results of refs. 
p5| , fl4f to four dimensions yields the appropriate weight 
function, which has a remarkably simple form in our co- 
ordinates, namely we minimize the following functional 
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where the wave-function is expressed in terms scaling 
functions (||), ^>k = ^i c i,i^k,i,i z - Due to our partic- 
ular choice of boundary conditions we consider only the 
terms for which I + l z is odd and l z — 3m, and truncate 
at I = l max = ka/2 + Al ka/2. 

We have computed three stretches of highly excited 
states. They consist of 7430, 1813 and 2362 consecutive 
eigenstates with 120 < ka < 235, 290 < ka < 300 and 
393 < ka < 400, respectively. The last two stretches 
comprise levels with sequential quantum numbers around 
20000 and 45000, respectively. The completeness of the 
series was checked by comparing the number of obtained 
eigenstates with the leading order prediction from the 
Weyl-formula d(k) = c(24 7 r 2 )~ 1 (afc) 2 , c w 0.51349. 

Fig. H shows that the nearest neighbor spacing distri- 
bution agrees very well with RMT predictions already for 
the lower energy spectral stretch 120 < ka < 235. The 
other series show well agreement, too. As for the long- 
range spectral correlations, the number variance S 2 (L) 
deviates from RMT predictions for interval length of 
more than ten mean level spacings which we believe is 
due to the parabolic-elliptic family of periodic orbits in 
the collinear manifold (Fig. |l|). The deviation from RMT 
decreases with increasing k. For the highest spectral 
stretch (ka ~ 400) the number variance increases linearly, 
Y?{L) fa T, 2 goc (L)+eL up to L < 250, with e ps 0.04. This 
finding is consistent with the model of a statistically inde- 
pendent fraction e of strongly scarred states flj4| . T, 2 (L) 
reaches its maximum and begins to oscillate at the sat- 
uration length L* which scales as L* oc k 2 in agreement 
with the prediction of ref. j26| . 

The length spectrum D{r), i.e. the cosine transform 
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of the oscillatory part of the spectral density d osc (k) = 
^ n i5(fc — k n ) — d(k), gives further information about 
long-range spectral fluctuations. For finite stretches 
of consecutive levels in the interval [ki , fe] one uses 
a Welsh window function w(k;ki,k2) = (&2 — k)(k — 
fci)/(6(fc2 — ^i) 3 ) m the- actual computation and ob- 
tains D(r) = J^ 2 dkw(k;ki,k 2 ) cos (kr) d osc (k) (see e.g. 

Fig. |] shows that orbits of length r = \/2a and 
its integer multiples cause dominant peaks in the length 
spectrum. 

To investigate the observed deviations from RMT pre- 
dictions in more detail it is useful to compute the in- 
verse participation ratio (IPR) of the wave functions in 
some basis pH| . In the case of billiards the use of the 
angular momentum basis (^|) is particularly convenient 
and suitable since periodic orbits correspond to sets of 

(n) 

isolated points within this representation. Let c\ ; de- 
note the expansion coefficients of nth eigenstate ^k n - 
We compute the IPR over a set of N consecutive eigen- 

states as IPR(U.) - ^E„ l<ffi I7(E„ \cQ | 2 ) 2 . The 
predicted RMT value for ideally quantum ergodic states 
is IPRgoc = 3. Fig. H shows the IPR for the two sets 
of eigenstates with 170 < k < 200 and 290 < k < 300, 
respectively. The agreement with RMT predictions is 
rather good in both cases. This confirms that the bil- 
liard under consideration is dominantly chaotic and er- 
godic. However, the IPR is slightly enhanced in the 
region around I = l z « fca/3. This is a robust phe- 
nomenon (present at all energy ranges), although the 
region of enhancement shrinks with increasing k. This 
finding is compatible with the expectation of uniform 
quantum ergodicity in the semi-classical limit. Note 
that the region I « l z corresponds to the vicinity of 
the collinear manifold. Note further that the orbits be- 
longing to the parabolic-elliptic family depicted in Fig. |l| 
have length V2a and angular momenta in the region 
l/ka — l z /ka G (1/2^6,1/^6)- This is precisely the 
region where the IPR exhibits its enhancement while the 
orbits' lengths coincide with the prominent peaks of the 
length spectrum in Fig. ||. 

Thus, the deviations from RMT predictions observed 
for the spectrum and for the wave functions are associ- 
ated with the family of parabolic-elliptic periodic orbits 
inside the collinear manifold. The special stability prop- 
erties of this family lead to scars in the wave functions 
of the quantum system. This is an exciting new type 
of scars of invariant manifolds and complements results 
previously found in a three-dimensional billiard |l4| and 
in interacting few-body systems [|l5| . Note that the fam- 
ily of parabolic bouncing ball orbits inside the collinear 
manifold does not cause statistically detectable scarring. 
The orbits of this family correspond to points in angu- 
lar momentum space with I = l z < kaj2\f% and do not 
exhibit an enhancement in the IPR since the classical 
motion is too unstable in their vicinity. 



In summary we have investigated an interacting three- 
body system realized as a billiard. Numerical results 
show that the classical dynamics is dominantly chaotic 
and no deviation from ergodic behavior is found. The 
spectral fluctuations of the quantum system agree well 
with random matrix theory predictions on energy scales 
of a few mean level spacings. However, wave function 
intensities and long ranged spectral fluctuations display 
deviations. These can be explained in terms of scars of a 
family of periodic orbits inside the collinear manifold. 
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FIG. 3. Integrated nearest neighbor spacing distribu- 
tion I(S) = J* S dsP(S) with the GOE value subtracted, 
I(S) — igoe(S), for the set of N = 7430 consecutive levels 
with 120 < k < 235 (full line). The dashed curve is the 
estimated statistical error a = y/l(S)(l-I(S))/N and the 
dotted curve is the difference for the commonly used Wigner 
surmise 7 w ig(<S') = 1 — exp(— 7rS 2 /4). 




FIG. 1. The motion inside the collinear manifold cor- 
responds to the motion inside a hexagonal billiard. The 
parabolic-elliptic family of periodic orbits is shaded. Three 
of its members and two limiting orbits are presented by full 
and dotted lines, respectively. 




FIG. 4. Length spectrum D(r) for the long spectral se- 
quence with 120 < k < 235. The dashed vertical lines denote 
the integer multiples of the period of the parabolic-elliptic 
family. 




FIG. 2. Billiard shape after separating off the center of 
mass motion and rotational degree of freedom. 
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FIG. 5. IPR(/,/ z = 3m), for 2052 states with 
170 < ka < 200 (left) and for 1813 states with 290 < ka < 300 
(right). Different levels of greyness are used for IPR values on 
consecutive intervals of width 0.5. The dominating light-grey 
corresponds to the interval [2.75, 3.25] around the RMT value. 
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